pip install pyLLEimport scipy
import plotly
import numpy
import matplotlib
import h5py
import prettytable
import matplotlib
import ipdb
If any errors show up due to a missing package please:
pip install <missing package name>using Pkg
Pkg.add("HDF5")
Pkg.update("HDF5")
Pkg.add("FFTW")
Pkg.update("FFTW")
Pkg.add("LinearAlgebra")
Pkg.update("LinearAlgebra")
Please note that if you are using this notebook, you need to make sure that you "trust" it to ensure the display of the figure. To do so, in the jupyter notebook, go to File → Trust this notebook. It should reload the notebook and the plotly figure will be able to be displayed now
In a python shell, spyder, a script, or in jupyter notebook, start by importing the package.
import pyLLE
We now define the resonator parameters. Here we will use a file TestDispersion.csv which is made of two columns: the first one is the azimuthal mode order (integer), the second is their corresponding frequency of resonance in Hz.
res = {'R': 23e-6, # ring radius in meter
'Qi': 1e6, # Intrinsic Q factor
'Qc': 1e6, # Coupled Q factor
'γ': 1.55, # Non-linear coefficient at the pump frequency
'dispfile': 'TestDispersion.csv', # frequency and corresponding azymuthal mode simulated previously
}
We now define the simulation parameters. Here we precise a linear detuning ramp of the pump from δω_init to δω_end relative to the pump mode angular frequency, mode closest to the defined pump frequency f_pmp. The simulation length Tscan is in unit of round trip, as it is more convenient in the Lugiato-Lefever formalism. It is important to notice that two parameters for the mode bandwidth have to be defined, μ_fit which determined the fit window of the raw data found in dispfile, and μ_sim which is the number of mode simulated in the LLE, hence could be larger than the fit mode through extrapolation
import numpy as np
sim = {'Pin': 150e-3, # Input power in Q
'Tscan': 1e6, # Length of the simulation in unit of round trip
'f_pmp': 191e12, # Pump Frequency
'δω_init': 2e9*2*np.pi, # Initial detuning of the pump in rad/s
'δω_end': -8e9*2*np.pi, # End detunin of the pump in rad/s
'μ_sim': [-74,170], # azimuthal mode to simulate on the left and right side of the pump
'μ_fit': [-71, 180], # azimuthal mode to fit the dispersion on the left and right side of the pump
}
In both the resonator and simulation dictionaries, the parameters can be called through their greek letters or through their equivalent latin names (e.g. μ → mu or δω → deltaomega). A translator dictionnary is implemented to translates every greek entries (see. self._greek)
Let's initialize the class
solver = pyLLE.LLEsolver(sim=sim,
res=res)
To plot and retrieve all the data of the dispersion, the method self.Analyze has to be called, resulting in the plot of the integrated dispersion Dint
dispfig = solver.Analyze(plot=True,plottype='all')
One can clearly see that because we simulate a larger window with the LLE (orange curve - LLE simulation) than the raw data (blue dots), we extrapolate outside of this region. One has to be carefull about this feature as ripple in the integrated dispersion can happen causing zero-crossings which are artefacts
A new attribute disp has been created which consists of a dictionary of the different value of the retrieve dispersion
solver.disp.keys()
One can solver the full temporal Lugiato Lefever equation :
$t_R \frac{\partial E(t, \tau)}{\partial t} = - \left(\frac{\alpha'}{2} - i\delta_0 \right)E + i \cdot \mathrm{FT}^{-1}\left[ -t_R D_{int}(\omega) \cdot \mathrm{FT}\left[E(t, \tau)\right]\right] + \gamma|E|^2 E + \sqrt{\theta}E_{in}$
The solver implemented in this package is based on a Julia core called through python. Hence to interface both language in an easy way, a hdf5 file is created in a temporary location with all the necessary data to solve the above LLE.
Hence, to solver the LLE here, we first need to setup the file:
solver.Setup()
And now we can solve the equation
solver.SolveTemporal()
It is important to note that a KeyboardInterupt will result in killing the Julia process
When done, we need to retrieve the data
solver.RetrieveData()
Here is a list of different method that ease the plotting of the quantity of interest, such as comb spectra and envelop of the electric field regarding the fast time. Although this is quite convenient for most of applications, one could use the raw complex electric field within the resonator. The different quantities saved by the solver are can be retrieved in the dictionary solver.sol
solver.sol.keys()
Everything in the solver.sol dictionarry is in SI units, except for the comb power, which is in dBm.
solver.sol['u_probe'], is the complex electric field relative to the fast time within the cavity, hence the complex number, probed by the solver. It is a matrix 1000 x [number of µ solved]solver.sol['Em_probe'] is the FFT of solver.sol['u_probe'], hence the complex spectra of the electric field. The comb spectra, within the cavity is given bu $10\times\mathrm{log10}(|\mathtt{Em\_probe}|^2)$ From Em_probe, it is easy to find the output comb power, given by Em_probe*sqrt(η), where η (in the Julia solver), is the coupling rate to the waveguide, given by Qcsolver.sol['Ewg'] corresponding to $\mathtt{solver.sol['Em_probe']}\times\sqrt{\eta}$, hence the complex spectra envelop of the electric field at the output of the coupled waveguide. solver.sol['detuning'] is the detuning ramp usedsolver.sol['ω'] is the angular frequencysolver.sol['ω'] is the optical frequencysolver.sol['theta'] is the position within the cavity, in radian (between -π and +π).It is import to note that all these quantity are subsampled in regards of the split-Fourier total step taken. Indeed, the number of step to solver the temporal LLE would be to large to efficienlty save or display the results. Hence, we subsample the steps taken by the solver at a fix number of 1000 (hence why the different matrices have a size of 1000x [number of µ solved])
We can plot the first figure that gives an overall idea of the behavior of the resonator
fig_sol = solver.PlotCombPower()